#!/bin/tcsh -f

# **********************************************************************
# TCtracker - Tropical Storm Detection
# Copyright (C) 1997-2008, 2021 Frederic Vitart, Joe Sirutis, Ming Zhao,
# Kyle Olivo, Keren Rosado and Seth Underwood
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 2
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
# 02110-1301, USA.
# **********************************************************************

# Useful directories
set LIBEXEC_DIR = @pkglibexecdir@
set PYTHONLIB_DIR = @pythondir@

# Default values
set year_beg = '0001'
set year_end = '0001'
set expname = 'Not Set'
set data_dir = /dev/null
set output_dir = /dev/null
set do_wsfc = y
set wind_th = 33.0
set matlab = 1
set ibtracs_dir = @ibtracs_dir@
set plots_only = n

if (${?IBTRACS_DIR}) then
  set ibtracs_dir = ${IBTRACS_DIR}
endif

set getopt_list = "b:e:n:O:i:o:pWw:Vh"
set argv = (`getopt $getopt_list $*`)
if ("$argv" == "--") then
   echo "To get help, use -h"; exit 1
endif
while ("$argv[1]" != "--")
   switch ($argv[1])
   case -b:
      set year_beg = "$argv[2]"; shift argv; breaksw
   case -e:
      set year_end = "$argv[2]"; shift argv; breaksw
   case -n:
      set expname = "$argv[2]"; shift argv; breaksw
   case -O:
      set ibtracs_dir = "$argv[2]"; shift argv; breaksw
   case -i:
      # Ensure data_dir is not a relative path
      set data_dir = `perl -e 'use File::Spec;print File::Spec->canonpath(shift)' "$argv[2]"`; shift argv; breaksw
   case -o:
      # Ensure output_dir is not a relative path
      set output_dir = `perl -e 'use File::Spec;print File::Spec->canonpath(shift)' "$argv[2]"`; shift argv; breaksw
   case -p:
      set plots_only = "y"; breaksw
   case -W:
      set do_wsfc = "n"; breaksw
   case -w:
      set wind_th = "$argv[2]"; shift argv; breaksw
   case -V:
      set print_version; breaksw
   case -h:
      set help; breaksw
   endsw
   shift argv
end
shift argv

if ($?print_version) then
   echo "tropical_storms_wsfc (@PACKAGE_NAME@) @VERSION@"
   exit 0
endif

if ($?help) then
   cat <<EOF
NAME
   tropical_storms_wsfc

SYNOPSIS
   tropical_storms_wsfc [-W] -b <year> -e <year> -n <experiment> -i <input directory> -o <output directory>

DESCRIPTION
   Creates figures and statistics from atmospheric data.

OPTIONS
   -W    Enables the use of w850 input data instead of wsfc.
   -b    First year to start processing data.
   -e    Last year to stop processing data.
   -n    Experiment name.
   -O    Location of ibtracs data.
   -i    Location of post-processed input files.
   -o    Location for storing the analysis data and figures.
   -p    Generate plots only
   -w    Set wind threshold
   -v    Print version information
   -h    Print usage information
EOF
   exit 1
endif

# Certain scripts/commands are GFDL specific
set gfdl = 0
if (`hostname -d` == "princeton.rdhpcs.noaa.gov") then
   set gfdl = 1
endif

if ($year_beg == "" || $year_end == "" || $expname == "" || \
    $data_dir == "" || $output_dir == "") then
   echo "Error: These options need values -- please consult the help screen (-h):"
   echo "Beginning year (-b): $year_beg"
   echo "Ending year (-e): $year_end"
   echo "Experiment (-n): $expname"
   echo "Input directory (-i): $data_dir"
   echo "Output directory (-o): $output_dir"
   exit 1
endif

# Check if data_dir exists
if ( ! -e $data_dir) then
   echo "Input directory '$data_dir' does not exist"
   exit 1
endif

set output_dir = $output_dir/atmos_${year_beg}_${year_end}/Zhao.TC

if ($gfdl) then
   source $MODULESHOME/init/csh
   module use -a /home/fms/local/modulefiles
   module purge
   module load gcp
   module load netcdf-c/4.7.3
   module load pyferret/7.6.3
endif

set tstorms_driver = $LIBEXEC_DIR/tstorms_driver@EXEEXT@
set hurricane_cat1to5_output = ${output_dir}/hurricanes_cat1-5_wsfc/
set hurricane_cat3to5_output = ${output_dir}/hurricanes_cat3-5_wsfc/
set hurricane_cat4to5_output = ${output_dir}/hurricanes_cat4-5_wsfc/

if ( "$do_wsfc" == "y" ) then
   set cyclones_dir = $output_dir/cyclones_wsfc
   set tstorms_dir = $output_dir/tropical_storms_wsfc/
   if ( ! -e $tstorms_dir  ) mkdir -p $tstorms_dir
else
   set cyclones_dir = $output_dir/cyclones_w850
endif

if ( ! -e $cyclones_dir  ) mkdir -p $cyclones_dir

set TSTORM_TEMP = $TMPDIR/$output_dir

if ($gfdl) then
   set TSTORM_TEMP = $FTMPDIR/$output_dir
endif
mkdir -p $TSTORM_TEMP
cd $TSTORM_TEMP

# For arbitrary years less than the smallest observation year (1980),
# we need to add the observation year to the current_year
# Discussed this with Ming.Zhao on 8/3/2013 - keo.
@ obs_first_year = ${year_beg}
@ obs_last_year = ${year_end}
if (${year_beg} < 100) then
   @ obs_first_year = ${year_beg} + 1980
   @ obs_last_year = ${year_end} + 1980
endif

@ last_completed_year = `find $cyclones_dir -name cyclones_\* -xtype f -exec basename {} \; | awk -F "." '{print $1}' | tail -1 | awk -F "_" '{print $NF}'`

if ($last_completed_year == 0) then
   @ next_year_to_process = $year_beg
else
   @ next_year_to_process += $last_completed_year + 1
endif

set file = `find ${data_dir} -name \*.nc -xtype f -printf "%f\n" | head -n 1`
set pp_type = `echo $file | awk -F "." '{print $1}'`
set file_first_year = `echo $file | awk -F "." '{print $2}' | awk -F "-" '{print $1}' | cut -c 1-4`
set file_last_year = `echo $file | awk -F "." '{print $2}' | awk -F "-" '{print $2}' | cut -c 1-4`
set start_month = `echo $file | awk -F "." '{print $2}' | awk -F "-" '{print $1}' | cut --complement -c 1-4`
set end_month = `echo $file | awk -F "." '{print $2}' | awk -F "-" '{print $2}' | cut --complement -c 1-4`

# If doing plots only, jump do plots
if ( $plots_only == "y" ) goto plots_only

# === Copy input files ===
echo "copying input files..."
# List of input variables
if ( "$do_wsfc" == "y" ) then
   set varlist = "u_ref v_ref"
else
   set varlist = "u850 v850"
endif
set varlist = "${varlist} vort850 tm slp"

# List of input files
cd ${data_dir}
set filelist = ()
foreach var (${varlist})
   #echo Processing variable ${var}
   foreach file (`ls -1 ${pp_type}.????010100-????123123.${var}.nc`)
      #echo Processing file ${file}
      # First year
      set yra = `echo $file | awk -F "." '{print $2}' | awk -F "-" '{print $1}' | cut -c 1-4`
      # Last year
      set yrb = `echo $file | awk -F "." '{print $2}' | awk -F "-" '{print $2}' | cut -c 1-4`

      # Does the time span of file overlap with desired period year_beg,year_end?
      # Note: use "bc -l" to remove any leading zeroes, otherwise numbers may be
      # interpreted as octals instead of decimals
      set n_yra = `echo "${yra}" | bc -l`
      set n_yrb = `echo "${yrb}" | bc -l`
      set n_year_beg = `echo "${year_beg}" | bc -l`
      set n_year_end = `echo "${year_end}" | bc -l`
      if ( ${n_yra} <= ${n_year_end} && ${n_yrb} >= ${n_year_beg} ) then
         set filelist = (${filelist} ${file})
      endif
   end
end

if ( $#filelist != 0 ) then
   mkdir -p $TSTORM_TEMP/tmp/
   if ($gfdl) then
      # Get files in chunks of 500 files
      echo ${filelist} | xargs -n 500 dmget
      set cp = gcp
   else
      set cp = cp
   endif
   ${cp} ${filelist} $TSTORM_TEMP/tmp/
else
   echo "Nothing to do, ending"
   exit 0
endif

echo "creating one year chunks, if needed..."
# Loop over years and create 1 year chunk files if necessary
cd $TSTORM_TEMP/tmp
foreach current_year (`seq $next_year_to_process 1 $year_end`)
   @ next_year = $current_year
   # Check if this is a non-January start date
   if ($start_month != "010100") then
      @ next_year = $current_year + 1
   endif

   # Padding year with initial zeroes, if needed
   set padded_current_year = `echo $current_year | awk '{printf "%04d\n", $0;}'`
   set padded_next_year = `echo $next_year | awk '{printf "%04d\n", $0;}'`

   set prefix = ${pp_type}.${padded_current_year}${start_month}-${padded_next_year}${end_month}

   foreach var (${varlist})
      test -e ${prefix}.${var}.nc
      set file_exist = $?
      if ( $file_exist ) then
         if ( $start_month == "010100" ) then
            # File with desired 1 year chunk length does not exist, create it.
            @ yra = $current_year
            @ yrb = $current_year + 1
            @NCRCAT@ -d time,"${yra}-01-01 06:00:00","${yrb}-01-01 00:00:00" \
                   -o ../${prefix}.${var}.nc ${pp_type}.*.${var}.nc
         else
            echo " ERROR: cannot create 1 year chunk files for non-January start date"
            exit
         endif
      else
         # File with 1 year chunk length exists, simply move it.
         mv ${prefix}.${var}.nc ..
      endif
   end
end

# Clean up tmp subdir
cd $TSTORM_TEMP
rm -rf tmp

# === Loop over years and process data ===
foreach current_year (`seq $next_year_to_process 1 $year_end`)
   echo "beginning process for year: $current_year"
   @ next_year = $current_year
   # Check if this is a non-January start date
   if ($start_month != "010100") then
      @ next_year = $current_year + 1
   endif

   # Padding year with initial zeroes, if needed
   set padded_current_year = `echo $current_year | awk '{printf "%04d\n", $0;}'`
   set padded_next_year = `echo $next_year | awk '{printf "%04d\n", $0;}'`

   set filename_prefix = ${pp_type}.${padded_current_year}${start_month}-${padded_next_year}${end_month}
   if ( "$do_wsfc" == "y" ) then
      set u_ref_filename        = ${filename_prefix}.u_ref.nc
      set v_ref_filename        = ${filename_prefix}.v_ref.nc
      set use_sfc_wnd = .true.
   else
      set u_ref_filename        = ${filename_prefix}.u850.nc
      set v_ref_filename        = ${filename_prefix}.v850.nc
      set use_sfc_wnd = .false.
   endif

   cat << eof > nml_input
 &nml_tstorms
   crit_vort  =  3.5E-5
   crit_twc   =  0.0
   crit_thick =  50.0
   crit_dist  =   4.0
   lat_bound_n =  70.0
   lat_bound_s = -70.0
   do_spline   = .false.
   do_thickness= .false.
 &end
 &input
   fn_u    = '$u_ref_filename'
   fn_v    = '$v_ref_filename'
   fn_vort = '${filename_prefix}.vort850.nc'
   fn_tm   = '${filename_prefix}.tm.nc'
   fn_slp  = '${filename_prefix}.slp.nc'
   use_sfc_wnd = $use_sfc_wnd
 &end
eof

   echo "beginning storm code..."
   cp $tstorms_driver tstorms_driver
   ./tstorms_driver < nml_input
   if ($? != 0) then
      echo ERROR: "tstorms_driver failed.  Giving up."
      exit 1
   endif
   if ($gfdl) then
      gcp nml_input $cyclones_dir
   else
      cp nml_input $cyclones_dir
   endif

   # To compare against obs data, we need to increase
   # the value of small start years (conversation with
   # Ming.Zhao 8/2/2013 - keo)
   @ obs_year = ${current_year}
   if (${current_year} < 100) then
      @ obs_year = ${current_year} + 1980
   endif

   echo "moving files to proper place..."
   if ( -e cyclones ) then
      mv cyclones cyclones_${obs_year}
      if ($gfdl) then
         gcp cyclones_${obs_year} ${cyclones_dir}
      else
         cp cyclones_${obs_year} ${cyclones_dir}
      endif
      echo "OUTPUT ADDED TO: " ${cyclones_dir}/cyclones_${obs_year}
   else
      echo "ERROR: JOB FAILED"
      exit
   endif
end

if ( "$do_wsfc" == "y" ) then
   echo "beginning hurricane trajectory code..."
   set trajectory_analysis  = $LIBEXEC_DIR/trajectory_analysis_csc@EXEEXT@

   cat << eof > nml_traj
 &input
   rcrit      =  1500.0
   wcrit      =  15.2
   vcrit      =  3.5e-5
   twc_crit   =  1.
   thick_crit =  0.
   nwcrit     =  3
   do_filt    = .true.
   nlat       = 70.
   slat       =-70.
   do_spline  = .false.
   do_thickness = .false.
   landmask_file = '$LIBEXEC_DIR/landsea.map'
   cmask_file = '$LIBEXEC_DIR/imask_2'
 &end
eof

   if ($gfdl) then
      gcp nml* $tstorms_dir
   else
      cp nml* $tstorms_dir
   endif

   set cyclone_files = `find $cyclones_dir -name cyclones_\* -xtype f`
   if ($#cyclone_files == 0) then
      echo "ERROR: No files to process."
      exit 1
   endif

   echo "running trajectory analysis for each cyclone file..."
   foreach cyclone_file ($cyclone_files)
      if (-e cyclones) then
         rm -f cyclones
      endif

      if ($gfdl) then
         gcp $cyclone_file cyclones
      else
         cp $cyclone_file cyclones
      endif
      set cyclone_filename = $cyclone_file:t
      set year = `echo $cyclone_filename | awk -F "_" '{print $NF}'` # also needed for plot scripts below
      $trajectory_analysis < nml_traj
      if ($? != 0) then
         echo "ERROR: trajectory analysis failed.  Giving up"
         exit 1
      endif
      mv ori_filt   ${tstorms_dir}/ori_${year}
      mv traj_filt  ${tstorms_dir}/traj_${year}
      mv trav_filt  ${tstorms_dir}/trav_${year}
      cat stats
      rm ori
      rm traj
      rm trav
      rm stats
   end

   plots_only:
   # Do plots
   # Determine the type of experiment [clim|obs]
   if ( ($obs_first_year < 1950 ) || ( $obs_first_year > 2005 ) ) then
      set experiment_type = clim
   else
      set experiment_type = obs
   endif

   echo "generating stats and plots"
   alias python3 @PYTHON@
   if ( $?PYTHONPATH ) then
      setenv PYTHONPATH ${PYTHONLIB_DIR}:${PYTHONPATH}
   else
      setenv PYTHONPATH ${PYTHONLIB_DIR}
   endif

   if ( $experiment_type == obs ) then
      # Snapshot plot, comparing observational data with model data
      python3 -m TCtracker.plots.snapshot -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   endif
   python3 -m TCtracker.plots.timeseries -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.by_latitude -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.by_longitude -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.by_region -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.seasonal_cycle -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.duration -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.jmaps_freq -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
   python3 -m TCtracker.plots.jmaps_ori -o $tstorms_dir $tstorms_dir $ibtracs_dir/obs $obs_first_year $obs_last_year $expname
endif
